Journal of Computer Graphics Techniques Vol. 8, No. 1, 2019 http://jcgt.org 


A Multiple-Scattering Microfacet Model for 
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Figure 1. Copper material (right) with and (left) without multiple scattering. 


Abstract 


This article introduces an extension to real-time image-based lighting models that incorpo- 
rates multiple scattering; it is suitable for real-time applications due to its low run-time cost. 
The main insight is that the precomputed integrals used for computing the single scattering 
term already contain all the information needed to simulate the remaining light bounces. The 
result is a technique that adds little overhead compared to its single scattering counterpart, 
but accomplishes perfect energy conservation and preservation. Even though the derivation is 
presented for a GGX bidirectional reflectance distribution function (BRDF), it can be trivially 
applied to other models as long as the split-sum approximation can be used to precompute the 


BRDF integral, or to find an analytical fit to it. 


1. Introduction 


The technique presented in this article builds on common image-based lighting mod- 
els for real-time applications. In Section 2, we introduce the basis of the technique by 
enforcing energy conservation on a perfect reflective surface. In Section 3, we gen- 
eralize the result to arbitrary conductors, taking Fresnel reflectance into account. We 
do so by reusing the precomputed data that was already used for the single scattering 
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term. Finally, in Section 4, we show how to extend the result to dielectric materials, 
and how the result both preserves energy at high roughness and prevents excess en- 
ergy at low roughness, yielding complete energy conservation and preservation. As 
shown in Figure 1, regaining this missing energy can have a strong visual impact on 
certain surfaces. 


1.1. Related Work 


Karis [2013] introduced the split-sum approximation as a method to precompute the 
single scattering BRDF integral for image-based lighting. His method uses prefiltered 
radiance with a precomputed integral for the BRDF response. He stores this integral 
as a two-component look-up table and uses it as a scale and bias for the Fresnel term. 

Kulla and Conty [2017] use the furnace test to compensate for the missing energy 
in the single scattering models and add an extra lobe to the BRDF based on that 
missing energy. The reasoning is analogous to the one presented in this article, but 
they compute several 2D and 3D look-up tables, and their results are intended for path 
tracing and not for real-time image-based lighting. 

Hill [2018a] refines the technique by Kulla and Conty and explores the expan- 
sion of the Fresnel term in multiple scattering events as a geometric series. Further 
[2018b], he uses a path-traced simulation to compute the energy emitted at every 
scattering event, which makes it more accurate than the model presented here, but 
leads him to precompute a look-up table for each event in the multiple scattering. His 
approach only explores analytical lights and doesn’t address the generalizations for 
image-based lighting. 


1.2. Furnace Test 


For energy to be conserved on a reflective surface, it must be true that for any incident 
light direction w;, the integral of reflected light along the viewing hemisphere equals 
the amount of incident light: 


[ BRDERr, Wi, Wo) COS Agdw, = 1, 
and given a reciprocal BRDF, the following must also hold: 
E(wo,r) = [ pepe. Wo, Wi) cos bidwi = 1. (1) 
Equation (1) defines the furnace test, which we’ll use several times in this article. 


The function E(wo,r) is also known as the directional albedo of a BRDF and is 
equivalent to integrating the BRDF against a uniform white environment. 
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Figure 2. Precomputed GGX first-bounce energy of a perfect reflector with horizontal axis 


cos(6,) and vertical axis roughness. 


2. Energy Preservation from Perfect Reflectors 


The simplest surface for which we can derive a multiple scattering model is the perfect 
reflector, which is a surface with reflectance F = 1. This eliminates the Fresnel term 
in the BRDF and simplifies integration. Integrating the single scattering BRDF (fss) 
against a solid white environment (Equation (2)), we obtain the energy that escapes 
the surface after the first bounce of light in a particular direction, the single scattering 
directional albedo Ess. Figure 2 shows the result of this integral using the same 
BRDF as Epic’s [2013]. In fact, it is equivalent to adding the scale and bias terms 
they propose: 


D( (Wi, Wot )G(wi, We, r) 


cos 6;dw;. (2) 


Ess (Wo,T 
i 4 cos 6; cos ĝo 


If energy is preserved, T 2 would be perfectly white, but since it only ac- 
counts for the first bounce of light, some areas are darker. By applying energy con- 
servation (Equations (3) and (4)), we can obtain the energy corresponding to the later 
bounces Ems: 

1 = Ess(wo, r) + Ems(wo,r), (3) 


Emal T] = f fms cos bidwi = 1 — Ess(wo, r). (4) 
Q 
We add an extra BRDF lobe to account for that missing energy: 
= f (fss + fms) cos 0; Lidwi = 1 fss cos 0; Lidwi + f fms cos ĝiLidwi. (5) 
Q Q Q 


On the right side of Equation (5), the first integral is the traditional single scattering 
term, which we can compute using Kulla’s split sum, and the second integral is the 
multiple scattering lobe. Similar to what we do for the single scattering term, we can 
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approximate this integral using a split sum (Equation (6)). The first term of the split 
is exactly Ems, and if we assume that energy participating in secondary scattering 
events is uniformly diffused, the second term can be approximated with the cosine- 
weighted irradiance, usually represented either as a low-resolution texture, or baked 
into spherical harmonics coefficients: 


Li 
f fms Li cos bidwi ~ f fms cos Qidwi f — cos 6;dw;. (6) 
Q Q Qa 


Notice that this approximation is actually exact for uniform environments, which 
explains the perfect furnace test in Figure 3. It’s still a good approximation for lights 
that cover a big portion of the hemisphere, where energy is more or less evenly dis- 
tributed, like sky domes, but not for narrow lights where most energy comes from a 
single direction (e.g., analytical point or directional lights). 

Equation (7) is the approximation for a perfect furnace test (image-based lighting) 
and Equation (8) provides the approximation for the split-sum calculation for perfect 
reflectors. 


i fms COS bids | — cos 6;dw; = (1 — Ea) f — cos ĝidwi, (7) 
Q QT QT 


Lo = Essradiance + (1 — E,,)irradiance. (8) 


Figure 3 shows a comparison between single scattering image-based lighting us- 
ing the split-sum approximation and the model with multiple scattering for perfect 
reflectors. The increased albedo is quite noticeable for high roughness values. 


Figure 3. Spheres with perfect reflector material (metallic=1, Fo = 1) with roughness from 0 
to 1; (above) an HDR environment and (below) a uniform environment (furnace test). The left 
half of the spheres uses the split-sum method, while the right half uses the multiple scattering 


approximation. 


3. Generic Metals, Fresnel Term 


For generic metals we need to take into account the Fresnel term. In order to reuse the 
split-sum approximation for perfect reflectors (Equation (8)), we'd like our equation 
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for uniform lighting to be in the form: 
E = FysEs5 + PmsEms:- (9) 


After the first bounce, light will be randomly scattered in all directions. We can 
model the distribution of light in these secondary bounces as having uniform energy 
in all directions, so it can be represented by an attenuated form of the cosine-weighted 
irradiance. Under that hypothesis, every secondary scattering event behaves the same 
way as the first one, except it takes this attenuated irradiance as its light source. This 
approximation makes sense because we already assume that incoming light is more or 
less evenly distributed when we use the split sum, and it will only get more diffused 
as it scatters by reflecting on randomly oriented microfacets. 

In much the same way as Hill [2018a], we can develop Fms Ems into a geometric 
series for the multiple bounces of light. On each light bounce, a fraction Eavg escapes 
the surface, and a fraction of 1 — Favg gets absorbed by the conductor, so only (1 - 
Braves) Fave is left for the next bounce, where Fy, is the cosine-weighted average of 
the Fresnel term. Notice we can’t approximate Fy, with Fss because while the latter 
depends on the viewing angle, the first represents attenuation for light coming from 
all directions, and so isn’t directional. However, Favg has an analytical solution when 
using the Schlick approximation (Equation (10) and (11)). 


a /2 6 
Fe 2r | (b= PA =a a =A) 
0 T 
1 
Favg = Fo + 5; (1 — Fo). (11) 
When Fo = 1, Equation (12) 
[0.6] 
E = Fzs Ess + > FzsEss(1 = Bigg) Eke (12) 


k=1 


and Equation (3) must be equivalent, so, 


(oe) (oe) 
1 = Eas + Y Bal — En)" = Eal — Eavg)*, 
k=1 k=0 


(oe) 
l= N Esl Ew = l 
2 ss( avg) 1 — (1 — Eag) 


Eavg = Ess. 
This implies that Favg is dependent on the viewing direction, as is Ess, meaning it 
takes more bounces for light to escape in some directions than others. By substitution 
of Eavg back into Equation (12), we have 
(1 — Eas) Fave 
> (1 ~~ Bias) Favg 


E = Fss Ess F Fss Ess I (13) 
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Figure 4. Precomputed GGX first bounce energy, separated into a term that depends on Fo 
(red) and an independent term (green). 


And combining Equations (9) and (13) we get 


(1 =T Ess)Favg 
1— (1 — Ess) Fave’ 


Fms Ems = Fms(1 _ Ess) = Fs Ess 


FssFavg Ess 
1— Fayg(1 — Ess) 


The problem with the previous equation is that we still don’t have an expres- 


Ema = (14) 


sion for Fss, but we can get one by integrating the single scattering BRDF. In Equa- 
tion (15), we do that using Schlick’s approximation to Fresnel: 


F; Ess =m fe fas (Wir Wo) (1 wo + wh)? cos 6;dw; 
F (wo, wh) 
e fesli Wo) (1 — wo- wh)” cos bidwi. (15) 
F (wo, wh) 


This is the integral that Epic Games computes into a look-up table for single 
scattering image-based lighting [2013], consisting of a scale (fa) and bias (fẹ) for the 
base reflectance Fo (see Figure 4), and we’ll refer to it as 


Fs Ess = Fofa + fo- (16) 


By substitution of Equation (16) into Equation (14), we get 


(Fo fa + fo)F: avg 


Fins = . (17) 
m1 = Fave(1 — Ess) 
So we can express the complete equation for metals as follows: 
F 
Lo = (Fofa + fo)radiance + (Fofa + fo) Favg Emsirradiance. (18) 


1— Fawg(1 — Ess) 
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// Common code for single and multiple scattering 


// vcoughness-—dependent Fresnel 
vec3 Fr = max(vec3(1.0 - roughness), FO) - FO; 


vec3 kS = FO + Fr * pow(1.0-ndv, 5.0); 


vec2 f_ab = textureLod(uEnvBRDF, vec2(ndv, roughness), 0).xy; 


vec3 FssEss = kS x f_ab.x + f_ab.y; 


float lodLevel = roughness * numEnvLevels; 


vec3 reflDir = reflect (-eye, normal); 


// Prefiltered radiance 

vec3 radiance = getRadiance(reflDir, lodLevel); 
// Cosine-weighted irradiance 

vec3 irradiance = getIrradiance (normal) ; 


// Multiple scattering 

float Ess = f_ab.x + f_ab-.y; 

float Ems = 1-Ess; 

vec3 Favg = FO + (1-F0)/21; 

vec3 Fms = FssEss*Favg/ (1-(1-Ess) *«Favg) ; 


// Composition 


return FssEss * radiance + Fms*Ems x» irradiance; 


Listing 1. GLSL code for generic conductors. 


Listing 1 shows an example implementation of image-based lighting for metals. 
Notice the overhead introduced by multiple scattering code is very low. Results can 
be compared in Figure 5. 


g 


<~ 


Figure 5. Spheres with copper material with roughness from 0 to 1; (above) in an HDR 
environment and (below) in a uniform environment (furnace test). The left half of the spheres 
uses the split-sum method for image-based lighting, while the right half uses the multiple 


scattering approximation. 
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4. Dielectrics 


With dielectrics there is a further term, corresponding to how the energy not reflected 
by Fresnel, diffuses under the surface and radiates back. Let Kg be the diffuse albedo 
of the surface, which models energy absorption under the surface, then 


E = Fss Ess + Fms Ems + KaEa. 


Eq is often computed as 1 — Fo. However, for grazing angles of incoming light, 
the Fresnel term will be higher (and less saturated), so less of that light should end up 
in the diffuse term. In order to get a better approximation of fq, we can observe that 
a perfectly white dielectric (Ką = 1) doesn’t absorb light, thus must radiate back as 
much energy as it receives: 


1= Fy5E'ss + FimsEms + Ea, 


Ea =1- Coe Ea Fms Ems). 


This fixes energy conservation for dielectrics, while still taking into account the 
multiple scattering specular lobe. It doesn’t explicitly model the effects of light scat- 
tering back and forth between the diffuse and specular interfaces, but since specular 
reflectance is usually quite low and unsaturated for dielectrics, radiated energy will 
eventually escape the surface unaltered, so the approximation holds. 

Figure 6 shows a furnace test on a grey background of the single-scattering (left) 
and multiple scattering (right) models for varying roughness. For high roughness, in- 
troducing multiple scattering recovers some missing energy, and for smooth surfaces, 
the excess energy towards the edges is eliminated. This is further illustrated in Fig- 
ure 7, where both models are compared to a reference photograph. Notice how the 
darker rim in the multiple scattering sphere (center), matches more closely the distri- 
bution of light in the reference. While the effect is certainly more subtle than with 
metals, removing the excess energy produces a more natural look. Listing 2 shows 
the code required to compute the proper Ką term for dielectrics. Again, the overhead 
over traditional diffuse IBL is quite small. 


Figure 6. Spheres with white dielectric material with roughness from 0 to 1 on a grey furnace 
environment. Left half of the spheres uses SS ibl with Eg = 1 — Fo, right half uses the MS 


approximation. 
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Figure 7. Close up of a smooth plastic ball. Single scattering model (left), multiple scattering 
model (center), and reference photograph (right). 


// GLSL code for dielectrics 
// Common code for single and multiple scattering 


// Roughness dependent fresnel 
vec3 Fr = max(vec3(1.0 - roughness), FO) - FO; 
vec3 kS = FO + Fr * pow(1.0-ndv, 5.0); 


vec2 f_ab = textureLod(uEnvBRDF, vec2(ndv, roughness), 0).xy; 


vec3 FssEss = kS x f_ab.x + f_ab.y; 


float lodLevel = roughness * numEnvLevels; 


vec3 reflDir = reflect (-eye, normal); 


// Prefiltered radiance 

vec3 radiance = getRadiance(reflDir, lodLevel); 
// Cosine-weighted irradiance 
vec3 irradiance = getIrradiance (normal); 


// Multiple scattering 

float Ess = f_ab.x + f_ab.y; 

float Ems = 1-Ess; 

FO + (1-F0)/21; 

vec3 Fms = FsskEss*Favg/ (1-(1-Ess) *Favg) ; 


vec3 Favg 


// Dielectrics 


vec3 Edss = 1 - (FssEss + Fms * Ems); 


vec3 kD = albedo » Edss; 


// Composition 


return FssEss * radiance + (Fms*Ems+kD) * irradiance; 


Listing 2. GLSLcode for generic dielectrics. 
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5. Conclusions 


Introducing the effect of multiple scattering back into image-based lighting can have a 
big impact on the appearance of materials. Aditionally, ensuring energy conservation 
and preservation makes materials more robust to changes in environment illumination, 
which is one of the benefits of physically based rendering in general. This article has 
shown that it is possible to achieve this effect with very little overhead, and that the 
visual impact can definitely be worth the cost. 
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